logo 
# Libraries
library(tidyverse)
library(sf)
library(USAboundaries)
library(sp)
library(leaflet)
library(plotly)
library(lubridate)
library(RColorBrewer)
library(knitr)
library(kableExtra)



Measurement threshold (1920 - 2019)

A minimum threshold of 10 distinct measurements with at least 1 measurement every 5 years was applied to filter for only regularly monitored wells








Active Management Areas

The 1980 Arizona Groundwater code established 8 Active management areas (AMA) in areas heavily reliant on mined groundwater. The goal of the AMAs is too achieve safe yield groundwater usage in these heavy groundwater pumping areas.




Joseph City AMA




Prescott AMA




Pheonix AMA




Harquahala AMA




Pinal AMA




Tucson AMA




Douglas AMA




Santa Cruz AMA










Buffer analysis on deep water wells

Investigation of wells below 300ft in which depth has continually decreased. The first plot shows all wells within a certain depth range that have had continuous measurement records for over 50 years. The following plots focus on a specific well in that range, and then a buffer was formed around these wells (12 - 20 km) to then see what is happening to other wells within the same area. The buffer analysis was used on individual wells that displayed dramatic and continual decreases in depth to water.





Range: 300 - 500 ft




Well 4573




Well 4740




Well 3472




Range: 400 - 500 ft




Well 4401




Well 4182




Range: 500 - 600 ft




Well 4993




Well 5730




Well 4953




Well 1853




Well 5743








Negative DTW measurements

A negative DTW reading may indicate groundwater emerging above the surface in the form of a spring. There are 9 wells that recorded negative depth to water readings. In the area of Chino Valley there are 6 wells with negative DTW to water readings all within a 5 km radius of one another. They other 3 wells are located in the far southeast corner of Arizona, nearing the US-Mexico border.

Negative Case Rules


  • Case 1A: Artesian to non-artesian (negative to positive sequentially)
    • maintained these values
  • Case 1B: Non-artesian to artesian (positive to negative sequentially, skeptical of these)
    • kept unless otherwise noted
  • Case 1C: Well with very low depth to water and fluctuation between slightly below ground surface and artesian
    • maintained these values
  • Case 1D: Well appears to have a pause in pumping activity (multiple years), returns to artesian, and then is pumped again
    • maintained these values
  • Case 2A: Measurements are sequential except for a random negative, therefore was likely a recording error
    • changed negative value to positive
  • Case 2B: Measurements are sequential except for periods of negative values that appear to be recorded incorrectly
    • changed negative value to positive
  • Case 2C: Measurements are sequential and within reason for well depth however all measurements for that well are negative and are likely recorded incorrectly
    • compared the well spatially with others to determine if artesian levels would be appropriate, and if not, changed to positive values
  • Case 3A: Measurements are sequential except for truly wonky number which doesn’t make sense
    • removed this observation
  • Case 3B: Measurements are sequential and then dataset changes to record 0 values only, or some other illogical number continuously, which does not make sense
    • removed these wells entirely





Case Count
1A 3
1B 0
1C 6
1D 0
2A 0
2B 0
2C 3
3A 0
3B 1

Negative DTW plots






Cases













Chino Valley












Functions

### PLOT AN INDIVIDUAL WELL BY WELL ID (DTW TIME SERIES)

### PLOT AN INDIVIDUAL WELL BY WELL ID (DTW TIME SERIES)
plotWell = function(df_time, num) {
  font = list(
    family = 'Courier',
    size = 15,
    color = 'white')
  label = list(
    bgcolor = '#232F34',
    bordercolor = 'transparent',
    font = font)
  well = df_time %>% filter(wellid == num)
  gg = ggplot(data = well, aes(x = date, y = dtw)) +
    geom_line(data = well, aes(y = dtw, col = wellid), size = 1) +
    ylim(max(well$dtw) + 100, 0) +
    labs(x = 'Year',
         y = 'DTW (ft)',
         col = 'Well') +  
    theme_bw() +
    theme(plot.title = element_text(face = 'bold',color = 'black', size = 18, hjust = 0.5), 
          axis.text.x = element_text(color="black", size=14), 
          axis.text.y = element_text(color="black", size=14), 
          axis.title.x = element_text(face="bold", color="black", size=16), 
          axis.title.y = element_text(face="bold", color="black", size=16), 
          panel.grid.major = element_line(colour = "#808080"),
          panel.grid.minor = element_line(colour = "#808080", size = 1))
  plot = ggplotly(gg, tooltip = c('x', 'y', 'wellid')) %>%
    style(hoverlabel = label) %>% 
    layout(font = font, 
           yaxis = list(fixedrange = TRUE)) %>% 
    config(displayModeBar = FALSE)
  plot
}

# PLOTS WELL DTW TIME SERIES WITHIN A RANGE (MIN, MAX)
plotRange = function(df_time, min, max) {
  font = list(
    family = 'Courier',
    size = 15,
    color = 'white')
  label = list(
    bgcolor = '#232F34',
    bordercolor = 'transparent',
    font = font)
  df_time = df_time %>% filter(dtw <= max, dtw >= min) %>% 
    arrange(desc(date))
  gg = ggplot(data = df_time, aes(x = date, y = dtw)) +
    geom_line(aes(y = dtw, col = wellid), size = 1) +
    ylim(max(df_time$dtw) + 50, min(df_time$dtw)) +
    labs(x = 'Year',
         y = 'DTW (ft)') + 
    theme_bw() +
    theme(plot.title = element_text(face = 'bold',color = 'black', size = 18, hjust = 0.5), 
          axis.text.x = element_text(color="black", size=14), 
          axis.text.y = element_text(color="black", size=14), 
          axis.title.x = element_text(face="bold", color="black", size=16), 
          axis.title.y = element_text(face="bold", color="black", size=16), 
          panel.grid.major = element_line(colour = "#808080"),
          panel.grid.minor = element_line(colour = "#808080", size = 1))
  plot = ggplotly(gg, tooltip = c('x', 'y', 'wellid')) %>%
    style(hoverlabel = label) %>% 
    layout(font = font, 
           yaxis = list(fixedrange = TRUE)) %>% 
    config(displayModeBar = FALSE)
  plot
  
}

# PLOTS WELL DTW TIME SERIES FROM DATAFRAME
plotMultipleWells = function(df_time) {
  font = list(
    family = 'Courier',
    size = 15,
    color = 'white')
   label = list(
    bgcolor = '#232F34',
    bordercolor = 'transparent',
    font = font)

  gg = ggplot(data = df_time, aes(x = date, y = dtw)) +
    geom_line(aes(y = dtw, col = wellid), size = 1) +
    ylim(max(df_time$dtw) + 50, 0) +
    labs(x = 'Year',
         y = 'DTW (ft)',
         col = 'Well') +  
    theme_bw() +
    theme(plot.title = element_text(face = 'bold',color = 'black', size = 18, hjust = 0.5), 
          axis.text.x = element_text(color="black", size=14), 
          axis.text.y = element_text(color="black", size=14), 
          axis.title.x = element_text(face="bold", color="black", size=16), 
          axis.title.y = element_text(face="bold", color="black", size=16), 
          panel.grid.major = element_line(colour = "#808080"),
          panel.grid.minor = element_line(colour = "#808080", size = 1))
  plot = ggplotly(gg, tooltip = c('x', 'y', 'wellid')) %>%
    style(hoverlabel = label) %>% 
    layout(font = font, 
           yaxis = list(fixedrange = TRUE)) %>% 
    config(displayModeBar = FALSE)
  plot
}

# IDENTICAL FUNCTION TO plotMultipleWells() BUT WILL PLOT ABOVE ZERO ON THE Y-AXIS
plotNegativeWells = function(df_time) {
  font = list(
    family = 'Courier',
    size = 15,
    color = 'white')
  label = list(
    bgcolor = '#232F34',
    bordercolor = 'transparent',
    font = font)
  gg = ggplot(data = df_time, aes(x = date, y = dtw)) +
    geom_line(aes(y = dtw, col = wellid), size = 1) +
    geom_hline(aes(yintercept = 0), size = 1) +
    scale_y_reverse() +
    labs(x = 'Year',
         y = 'DTW (ft)',
         col = 'Well') + 
    theme_bw() +
    theme(plot.title = element_text(face = 'bold',color = 'black', size = 18, hjust = 0.5), 
          axis.text.x = element_text(color="black", size=14), 
          axis.text.y = element_text(color="black", size=14), 
          axis.title.x = element_text(face="bold", color="black", size=16), 
          axis.title.y = element_text(face="bold", color="black", size=16), 
          panel.grid.major = element_line(colour = "#808080"),
          panel.grid.minor = element_line(colour = "#808080", size = 1),
          legend.title = element_text(colour="black", size=16, face="bold"),
          legend.text = element_text(colour="black", size=10, face="bold"))
  plot = ggplotly(gg, tooltip = c('x', 'y', 'wellid')) %>%
    style(hoverlabel = label) %>% 
    layout(font = font, 
           yaxis = list(fixedrange = TRUE)) %>% 
    config(displayModeBar = FALSE)
  plot
}

# RETURNS A TIME SERIES PLOT OF WELLS WITHIN A BUFFER DISTANCE (km)
plotBuffer = function(df, id, buffer) {
  font = list(
    family = 'Courier',
    size = 15,
    color = 'white')
  label = list(
    bgcolor = '#232F34',
    bordercolor = 'transparent',
    font = font)
  well = df %>% filter(wellid == !!id)
  buff = st_buffer(df[well, ], buffer)
  nearby =st_intersection(df, buff)         
  df_time = join_time %>% filter(time_span > 50, wellid %in% nearby$wellid) 
  
  gg = ggplot(data = df_time, aes(x = date, y = dtw)) +
    geom_line(aes(y = dtw, col = wellid), size = 1) +
    ylim(1000, 0) +
    labs(x = 'Year',
         y = 'DTW (ft)',
         col = 'Well') +
    theme_bw() +
    theme(plot.title = element_text(face = 'bold',color = 'black', size = 18, hjust = 0.5),
          axis.text.x = element_text(color="black", size=14),
          axis.text.y = element_text(color="black", size=14),
          axis.title.x = element_text(face="bold", color="black", size=16),
          axis.title.y = element_text(face="bold", color="black", size=16),
          panel.grid.major = element_line(colour = "#808080"),
          panel.grid.minor = element_line(colour = "#808080", size = 1))
  plot = ggplotly(gg, tooltip = c('x', 'y', 'wellid')) %>%
    style(hoverlabel = label) %>% 
    layout(font = font, 
           yaxis = list(fixedrange = TRUE)) %>% 
    config(displayModeBar = FALSE)
  plot
}


# PLOTS A BUFFER AROUND SPECIFIED WELL AND LOCATES OTHER WELLS INSIDE THE BUFFER AREA
buffer_fun = function(df, well, buff, state) {
  buffer = st_buffer(df[well,], buff)
  near1 = st_intersection(df[,], buffer) %>% filter(measurement_dist >= 10)
  
  plot = ggplot() + 
    geom_sf(data = state) +
    geom_sf(data = buffer, fill = NA) + 
    geom_sf(data = near1, col = "red", size = .5) + 
    labs(caption = paste(nrow(near1), 'wells')) +
    theme_void() +
    theme(plot.caption = element_text(size = 22, face = "bold", hjust = 0.5))
  print(plot)
  return(near1)
}

# RETURNS DATAFRAME OF WELLS IN DTW RANGE (MIN - MAX)
dtw_range = function(df, min, max) {
  df = df %>% filter(dtw <= max, dtw >= min) %>% 
    mutate(sd = sd(dtw)) %>% 
    arrange(desc(date))
  return(df)
}